function points=ElixhauserPointSystemLogit(elix,y,alpha)


%%%%%
% This aggregates the 30 Elixhauser comorbidities into a point for a particular choice of binary outcome (e.g., in-hospital mortality).
% 
% Input:
% elix - Nx30 matrix containing binary Elixhauser comorbidities, where N is the number of subjects.
% y - column vector containing the binary outcome.
% alpha - statistical significance level that determines which comorbidities are included in point calculation
%
% Output:
% points - column vector containing points for the corresponding rows in elix.
%
% Written by Joon Lee, 2011
%%%%%


idxincl=1:size(elix,2);
idxincl(all(elix==0) | all(elix==1))=[];
finish=false;

while ~finish && ~isempty(idxincl)
    [b,dev,stats]=glmfit(elix(:,idxincl),y,'binomial','link','logit');    
    p=stats.p(2:end);
    idxincl(isnan(p))=[];
    b([false; isnan(p)])=[];
    p(isnan(p))=[];    
    if ~all(p<alpha)
        [C,I]=max(p);
        idxincl(I)=[];        
    else
        finish=true;
    end    
end

points=zeros(size(elix,2),1);
if ~isempty(idxincl)
    points(idxincl)=round(b(2:end)./min(abs(b(2:end))));
end
